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Abstract. We construct a family of root-finding algorithms which combine knowledge of the 
branched covering structure of a polynomial with a path-lifting algorithm for finding individual roots. 
In particular, the family includes an algorithm that computes an e-factorization of a polynomial of 
degree d which has an arithmetic complexity of O(d(logfi)^| loge| -I- (P{\ogd)^). At the present time, 
this complexity is the best known in terms of the degree. 
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Introduction. The problem of devising optimal methods for numerically approx- 
imating the roots of a polynomial has been of interest for several centuries, and is far 
from solved. There are numerous recent works on root-finding algorithms and their 
cost, for example, the work of Jenkins and Traub ||JT70|| , Renegar [[Ren87|| , Schonhage 
Sch82|| , and Shub and Smale |[SS85| , [SS86| , |Sma85|| . This list is far from complete; the 



reader should refer to the aforementioned papers as well as [ pH69|| for more detailed 
background. The work in this paper is most closely related to that of Smale. 

Our algorithm computes an approximate factorization of a given polynomial (that 
is, it approximates all the roots). In constructing it, we combine global topological in- 
formation about polynomials (namely, that they act as branched covers of the Riemann 
Sphere) with a path-lifting method for finding individual roots. Utilizing this global 
information enables us to use fewer operations than applying the path-lifting method 
to each root sequentially. 

Renegar's algorithm in |[Ren87|| approximates all d roots of a given polynomial 
using logfi -|- (i^(log (i)(log I loge|)^ arithmetic operations in the worst case. He 
has shown that the factor of log|loge| in the complexity is the best possible if one 
restricts to the operations -|-, — , x, and -i-. This algorithm has a component (the 
Shur-Cohn algorithm) which requires exact computation and so is not suitable for an 
analysis of bit complexity, that is, one which accounts for rounding errors introduced 
by finite precision. In |Pan87|| , Pan gives a different algorithm which slightly improves 
the complexity to 0(^(i^log(ilog | loge|^; this algorithm also operates effectively as a 
parallel algorithm. 



Schonhage |Sch82| gives an algorithm which produces an e-factorization with a 
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bit complexity of 0(^(1^ \ogd + (f\ loge|^ log((i| loge|) loglog((i| loge|), via the "splitting 
circle method". Note that the customary parameter for bit length of the coefficients 
does not appear in the complexity. This is because, as Schonhage states, for fixed 
degree d and output precision e, there is a number sq for which "the input [coefficients] 
a^, can be restricted to complex integer multiples of 2"*'^ without loss of generality." 



In |[Ren87|| , it is stated that Schonhage believes that, if exact arithmetic is used, this 
method "should yield a complexity bound [in e] of log | loge|), most likely with 
a < 3." 

Smale's path lifting algorithm in |pma85|] approximates a single root of the polyno- 
mial with a worst case arithmetic complexity of 0(^d{\ogd)\ loge| + (i^(log(i)^^, and an 

average complexity of 0(^d^ + d\ loge|^. One good feature of this line of work is that 
it is stable under round-off error. In |[Kim89a[ , Kim shows that if / and /' are com- 
puted with relative error 10-3 until an approximate zero (see Section |1.2|) is reached, 
then the algorithm behaves exactly the same. A recent series of papers by Shub and 
Smale ||SS934 ^S93b| , ^S93c| , |SS93d|| generalizes the path lifting algorithm to systems 
of homogeneous polynomials in several variables. 

The algorithm presented here exploits the branched covering structure of a polyno- 
mial to choose good starting points for a variant of Smale's algorithm, and we obtain 
a worst case arithmetic complexity of 0(^(i(log(i)| loge| + (i^(log(i)^^ to compute an 
e-factorization. In a subsequent paper, we shall compute the bit complexity of this al- 
gorithm. Because of the stability mentioned in the previous paragraph and the ability 
to exploit bounds on the variation of / and /', we hope to achieve results comparable 
to Schonhage's. 

At first glance, it may appear that our complexity results are inferior to some of 
those above in terms of e. However, in practice there is usually a relationship between 
the degree d and the desired precision e; if we have e > 2~'^, then the complexity of 
our algorithm compares favorably with all of those mentioned above. Furthermore, our 
algorithm is quite simple to implement and is numerically very stable. 

Our algorithm is suitable for some amount of parallelization, but has a sequential 
component of 0{d+ \ loge|) operations. However, we think of this algorithm as acting 
on d points simultaneously, and techniques which evaluate a polynomial at d points (see 
BM75 |, for example) are used to cut the cost involved. Of course, the algorithm can be 
implemented on a sequential machine while still taking advantage of these techniques. 
In fact, evaluation of the polynomial is the only point at which we at which we need to 
use asymptotic estimates to achieve the stated complexity; the other places where we 
use asymptotic estimates are only for ease of expositon. 

The reader should also see the papers ||BFKT88| , pT90| , [Nef9CI|| for fully parallel 
algorithms for solving polynomials with integer coefficients. In ||BFKT88|] , it is shown 
that if all roots of the polynomial are real, this problem is in NC. Neff extends this 
result to allow complex roots in ||Nef9CI| . 

This paper is structured as follows: In Chapter [l|, after some background material, 
we recall the "path lifting method" of Smale and present a version of the relevant 
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theorem (our Theorem |1.5| ) which improves the constants involved somewhat. We then 
discuss how we can exploit the branched-covering structure of a polynomial to choose 
initial points for the algorithm, many of which will converge to roots. We close the 
chapter with a brief explanation of how to construct families of algorithms which locate 
d/n roots at a time, for various values of n. 

Chapter |^ presents an explicit algorithm for a specific family, which locates d/2 
points at a time. Our main theorem. Theorem states that this algorithm always 
terminates with an e-factorization of the input polynomial, and gives a bound on the 
number of arithmetic operations required in the worst case. As a corollary, the algorithm 
can be used to locate all d roots of the polynomial to within e with a complexity of 
0(^(i^(log(i)| loge| + (i^(log(i)^^. In the subsequent sections, each component of the 
algorithm is analyzed, and the relevant lemmas are proven. Finally, we tie all the 
components together and prove the main theorem. 



1. Preliminaries. 

1.1. Root and coefficient bounds. Given a polynomial (j){z) = I]f=o'^j-^*5 with 
Oi G C, it is our goal to determine an approximate factorization of 0, that is, approxi- 
mations to the actual roots of so that \\(f>{z) — Yl{z — ^t) \\ < e. The norm we shall 
use here is the max-norm, that is, ||0|| = max \ aj\. A related problem is to ensure that 
16 ~ 61 < f'; there are well-known estimates giving the relationship between e and e', 
so solving one problem essentially solves the other. 

In order to have an estimate on the complexity of a root-finding algorithm, we 
need a compactness condition on the space of polynomials. This can be done either by 
placing conditions of the location of the roots or on the coefficients; such bounds are 
interrelated. 

Since our goal is to minimize a functional norm, it seems most natural to place our 
conditions on the coefficients. We shall assume our input polynomial is an element 
of the family 

p^(l) = {^'^ + ^ ajz\ with \aj\ < 1}. 

j=0 



This is the same space as used by Smale and others [[SS85| , |Sma81| , [Sma85|| . One can 
always transform an arbitrary polynomial into an element of ^^(1): if p{z) = Y,i=o 
and B = max \bj/bd\^^^'^-^\ then p{Bz)/B^ e Vd{l). 

One should not confuse this family with the degree d polynomials whose roots are 
in the unit disk, although unfortunately this space is also often denoted by Vd{l) (for 
example, in |[I'ri9(]|| and ||Ren87|| ). 



There are a number of estimates which relate the coefficients of a polynomial to 
a bound on the modulus of the zeros (see ||Hen74| ] or ||Mar66|| , for example). Such an 



estimate is important to us, since although membership in 'Pd(l) is not preserved under 
deflation (division of factors), bounds on the modulus of the roots are. We state one 
such bound here (Corollary 6.4k of [ IIen74|] ) : 
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Lemma 1.1. All the zeros of the polynomial z'^ + ^2"^^^ a jZ^ lie within the open disk 
with center and radius 

2 max \aA^'^^-^\ 

0<j<d 



As an immediate consequence, we see that the roots of a polynomial mVd{l) lie within 

»2. 



1.2. Approximate zeros. Our algorithm uses a path lifting method (see below) 
to get close to the roots of our polynomial, and then uses the standard Newton's method 
to further refine these approximations. This is done because Newton's method converges 
very quickly in a neighborhood of a simple root, but can fail for some initial points 



outside this neighborhood. One of the authors ||Sut89|| has shown how one can guarantee 
convergence of Newton's method, but a bound on the arithmetic complexity has not 
been computed. Instead, we use the more certain path lifting method as described in 
Section |1.3| ; this allows an explicit computation of the complexity. 

Following Smale |^ma81|| , we call a point zq an approximate zero if Newton's method 
converges rapidly (that is, quadratically) when started from zq. Such terminology is 
reasonable, because given such a point, we can quickly obtain an approximation of a 
root to arbitrary precision. 

Definition 1.2. Let f be a polynomial and let Zn be the ra*^ iterate under Newton's 
method of the point zq, that is, Zn = Zn-i — f{zn-i)/f'{zn^i). Then we say that zq is 
an approximate zero of f if, for all n > we have 



2n 



ko - CI, 



for some root ( of f . 

Notice that this definition is never satisfied in the neighborhood of a multiple root 
of /, since the convergence of Newton's method is asymptotically linear there. In our 
algorithm, we perturb the polynomial slightly to ensure that we always have simple 
zeros. Refer to Section and Section for more details. 

Kim [ [Kim88|| and Smale ||Sma86|| have developed readily tested criteria for deter- 
mining, based on the values of the derivatives f'^^\z), when a point z is an approximate 
zero. These can be extended to a much more general setting, namely for / a mapping 
between Banach spaces. The following is essentially Theorem A of [^ma86|| : 

Lemma 1.3. Let 



z) = max 

k>l 



f{z) 



k\f'{z) 



i/(fc-i) 



If af{z) < \, then z is an approximate zero of f . 
We will find the following also very useful. 
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Lemma 1.4. Let f{z) be a polynomial and z be a complex number so that f'{z) ^ 0, 
and let Rf{z) be the radius of convergence of the branch of the inverse f~^ which takes 
f{z) to z. If 

Rf{z) 10' 

then z is an approximate zero of f . Furthermore, if we have \ f{z)\/ Rf{z) < 1/32, then 
«/(^) < 1/8. 

Remark. If Smale's mean value conjecture holds (see |^ma81| , |Tis89|| ), then the 
hypotheses of the lemma imply that otf{z) < 1/8. 



Proof. The first result is a consequence of the proof of Theorem 4.4 of ||Kim88 
and the second is an immediate consequence of Corollary 4.3 of the same paper. This, 
in turn, uses the Extended Loewner's Theorem in |[Sma81|| . □ 

1.3. The path lifting method. Here we review the path lifting method, which 
forms the core of our iteration scheme. This method is sometimes referred to as a 
"generalized Euler method" or "modified Newton's method"; we prefer the term "path 
lifting method" as it is the most descriptive (to us, anyway). This method appears in the 
work of Steven Smale |^ma85|| , although the version we present here is slightly different 



and we present another proof of the relevant theorem, which is quite simple. It should 
be emphasized that the path lifting method, like Newton iteration, is an algorithm for 
finding a single root of a polynomial; we discuss how to combine these to find all roots 
in Section |l]5] below. 

We think of a polynomial / as a map from the source space to the target space] 
that is, / : <Csource ^target- Givcu an initial value zq in the source space, we connect 
its image wq = f{zo) to in the target space, and then lift this ray under the proper 
branch of f~^ to a path connecting zq with a root ( of /. Of course, we don't explicitly 
know this inverse, but if the path in the target space stays well away from the critical 
values of /, the local inverse map f~^ is well-defined on a neighborhood of the ray. Even 
if the path does contain critical values, a local inverse can still be defined for some zo- 
See Section |L 



z„ 





source space 

Fig. 1.1. The source and target spaces in the path lifting method. In the source space, each Zi is 
indicated by a black dot, and f{zi) is indicated by a black dot in the target space. Similarly, the Wi are 
indicated by tick marks in the target space, and f~^(wi) by ticks in the source space. 
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The basic idea of the path hfting method is to take a sequence of points Wn along 
the ray in the target space, with wq = f{zo). We then construct a sequence of points 
Zn in the source space so that f{zn) is near w„ in the target. This is done using a single 
step of Newton's method to solve f{z) = Wn with Zn-i as the starting point. That is, 

_ _ f{Zn-l) - Wn-l 

f'{Zn-l) • 

This construction will converge to a root ( in the source space if there is a wedge 
about the ray in the target space on which there is a well-defined branch of the inverse 
/^^, and if the Wn are chosen properly (in a way which depends on the angle of the 
wedge) . The larger the wedge about the ray, the faster the method converges. We now 



state the exact theorem, although we shall defer the proof until Section |2T3 . 
Notation. By a wedge Wa,w, we mean the set 

{z < \z\ < 2\w\, arg w — A < arg z < argw + A}. 

Theorem 1.5. Suppose that the branch of the inverse is analytic on a wedge 
yVA,wo! ''^ith < A < n/2, and let h < ^^f^- Suppose also that |/(zo) — wo\ < h\wo\/2, 
and define 

Wn= [l- h) Wo, Zn+1 = Zn 7 . 

Then |/(z„) - m„| < h\wn\/2 and Zn+i G /^^^ (>Va,«;„) • 

It should be noted that this theorem is a slight improvement of Smale's Theorem A 
|Sma85|| . His proof is valid for all angles, but is stated only for A = 7c/12 with 



m 



h = 1/98, and for A = n/A with h = 1/32. For A = 7r/12 we can take take h = 1/74, 
and for A = it/A, h = 1/27 is adequate. 



1.4. Branched covers, inverse functions, and all that. If / is a polynomial 
of degree d, then / : C — >■ C is a branched covering with branch points at the critical 
points 6i of /. If 2; is a regular point of /, that is, f'{z) 7^ 0, then there is a well-defined 
inverse function so that (/(-z)) = z. 

In any neighborhood of a critical point of /, there cannot be a single valued inverse; 
however, the behavior at such points is well understood. Let ^ be a critical point of 
multiplicity k — 1. Then we have 

f(z)-fie) = iz-efgiz), where 9(6)^0. 

One can then define k branches of the inverse which are analytic on a small slit disk 
about f{9). We may, of course, choose any slit which connects f{9) to the boundary 
of the disk. The reader is referred to a complex analysis text for further details (for 
example, see ||Ahl79|| ). 
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Let {Q} be the d roots of /, represented with multiphcity. If Q is a simple root, 
denote by /^^ or f^^ the unique branch of the inverse of / which takes to Q. On 
the other hand, if Q is a multiple root, we let /^T^ = fj'^ be one of the branches of 
the inverse that take to Q, taking care to account for all such branches exactly once. 
This can be done, since if C is a root of multiphcity A; > 2, it is also a critical point of 
multiplicity k — 1, and so there are k branches of the inverse. 

We now analytically continue each of the to a maximal starlike domain Q,j in 
the target space; that is, we attempt to extend each along open rays from 0. When 
doing this, it is useful to think of the target space as consisting of d copies of C, with a 
single fj'^ associated to each one. When does the analytic continuation fail? Precisely 
when the inverse image of a ray encounters a critical point of /. Refer to Figure 1.2. 




Fig. 1.2. The flj in the target space (on the right), and the corresponding source space, for a 
degree 5 polynomial with a double root (black dot inside a white dot) and a critical point of multiplicity 
2 (double white dot). The other critical point is marked by a single white dot, and the roots by black 
dots. The cuts in the flj are represented by black radial lines. 

At this point, it may be useful to consider the Newton vector field given by 

z^-f{z)/f{z). 

Let (fit{z) be the solution curve with initial condition (fio{z) — z. Then we have 

fMz))=e-'f{zy, 

that is, / maps solution curves of the Newton vector field to rays in the target space. 
Notice that the singularities of the vector field occur precisely at the points where 
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f'{z) = 0. The solution curves (ft{z) which have singularities play an important role 
here: they divide the source space into regions on which / is injective. Refer to |^hu83 



|STW88| , |Sma85|| for more details on the behavior of the solution curves. Applying the 
path lifting method can be viewed as attempting to follow the solution curves to the 
flow (pt- 

When constructing the we continue radially outward until a critical point 
6 is encountered in the source space. We then exclude the ray {rf{6) | r > 1} from Qj, 
and continue by moving along rays which avoid the cut. Notice that when we encounter 
a critical point 6 of multiplicity k — 1, we need to slit at most k of the flj starting at 
f{0). Also, note that some of the Qj may already be slit at f{0), since there may be 
another critical point whose image lies on the same ray. 

We now count the number of such cuts: / has d—1 critical points (with multiplicity), 
and a critical point of multiplicity k — 1 can cause at most k cuts. This means we have 
at most 2{d — 1) cuts, distributed through the d copies of the target space. Note that 
if flj contains some wedge W, then is analytic on W. The following counts the 
number of Qj which contain wedges of a given size. 

Lemma 1.6. Let m be an integer, and divide C into m wedges 



2mT 



< argw < 



2(n + l)TC 



m 



m 



} 



= 0, . . . , m — 1. 



For each wedge Wn,m, ^et A^(n, m) he the number of Vtj which contain the sector Wn^m, 
and let N{m) = max N{n,m). Then 

0<n<m 



N{m) > d 



2{d-l] 



m 



Proof. Since we have 2{d — 1) cuts and m wedges Wn, at least one of the wedges 
has no more than ^iitill cuts. Since there are d f2,s, we have the result. □ 

Corollary 1.7. N{d) = d, N{3) >-, and N{4:) > -. 

3 2 

Proof. Application of the formula above gives the result for A^(3) and A^(4), and 
yields N{d) > d — 1. However, since each critical point causes at least two cuts, it is 
not possible to have a wedge cut only once. Thus, the value d — 1 is not permissible for 
N{d), giving N{d) = d. □ 



1.5. Families of root-finding metliods. We now use the result of Lemma |L6 
to construct families of root-finding algorithms. Recall that the the modified Newton 
method described in Section |1.3| works when there is an Qj containing a wedge about 



our initial value f{zo); the larger the wedge, the faster the method converges. 

For each family, we start with md points in the source space placed around a circle 
which contains all the roots. We think of this as m sets of d initial points, and choose 
them so that the image of each set lies well inside each of the m sectors Wn^m in the 



target space. Then by Lemma |1.6|, one of the m sets of points will contain at least 
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N{m) elements whose images are each in a "good wedge", that is, they he in some flj. 
As a consequence, iterating these points under the path hfting method will locate at 
least N{m) roots of the polynomial. 

Particular families of interest are m = d, which gives the algorithm discussed in 



Kim89b|| , and m = 4, on which we focus our attention here. The basic idea of all 
of the algorithms is this: obtain md "good" initial points and apply the path lifting 
method to d of them at a time. If, after a prescribed number of iterations, we have found 
approximation to at least N{m) roots (counting multiplicity), we deflate the polynomial 
(that is, divide out the approximated roots) and repeat the process on the result. If 
not, we try again with the next set of d points. Note that we are guaranteed success 
by the time we try the m*'^ set. 

The remainder of the paper consists of a detailed description and analysis of the 
algorithm for m = 4. Most of what follows can be readily adapted to the other families 
as well. 

2. A root-finding algorithm. 

2.1. Statement of the algorithm and main theorem. Here we present our 
root-finding algorithm for the family m = 4. The presentation is structured as a main 
routine and several subroutines, which do most of the work. 

Notation. Throughout this chapter, we shall denote matrices, vectors, and sets 
in uppercase calligraphic type, and their elements in subscripted lowercase type. For 
example, Xj is the j^^ element of the vector X. We shall also use the notation \_x\ to 
denote the least integer in x, sometimes also called f loor(a;). 

The main routine merely inputs the desired polynomial and precision, rescales it so 
the roots lie in the disk of radius 1/2, then repeatedly calls a subroutine to halve the 
number of unknown roots (counted with multiplicity) and deflate. We do the rescaling 
in order to easily bound the error introduced by the FFT deflation. The set A contains 
all the approximations found by the i^^ stage. 

Note that the algorithm is given for an arbitrary monic polynomial, since only 
minor changes are required to normalize the input polynomial. If it is assumed that the 
input polynomial is already in Pd(l), we can take fo{z) = 0(42;)/4'^ and r = 32e/7'^+^. 

Main Routine 

Input monic polynomial 4>{z) = J2i=o aiz'' and desired precision e. 
Let /o(z) = (l){Kz)/K'^, with K = 'i max|aj|d^. 



Let r = ^ {Anr\ 
Let i = and A = 
While #(A) < deg((/>) 

Let (/i+i,Aj) = get-half -roots-and-def late(/j, r). 

Let A = A U Aj. 

Increment i. 
End While 
Output Kk. 
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The function get-half -roots-and-def late takes as input a normahzed polyno- 
mial / and precision r. It returns a set of points yj which approximate at least half of 
the roots of / (with multiplicity) and a new polynomial / which we obtain by deflation. 
These satisfy \\f{z) — f{z)Y{{z — yj)\\ < 2r. We should point out here that we are 
actually finding approximate zeros of / — r , where |r| = r, which depends on e. When 
the translation is in the proper quadrant, this will ensure that the relevant roots of 
/ — r are simple, so that we have approximate zeros in a neighborhood. This allows 
us to obtain the right number of approximations to a multiple root, without worrying 
about winding number arguments or the like. We emphasize again that r is chosen as 
a function of e, and is small enough that the approximation polynomial has negligible 
errors in the non-constant terms. 

The matrix Z consists of 4 rows of d "good" initial conditions, with \zj^k\ = 3/2 
and aigf{zj^k) ^ 27rij'/4. We use Zj to represent the j^^ row. 

function get-half -roots-and-def late(/, r) 

Let d = deg /. 

Let Z =choose-4d-good-initial-points(/). 
For j = 1 to 4 do 

Let 3^ = iterate-PLM(/, Zj, j, r). 

Let = / - re^^^J'/^. 

Let X = select-approx-zeros(^, 3^). 

Let W = polish-roots('(/', A', r). 

Let V = weed-out-duplicates('(/', W). 

If #(V) > d/2 then 

Let / = deflate(V', V). 
Return (/, V) 

End if 
End for 

Print "We have proven this statement will not be reached. Check your code." 
Abort. 

The function choose-4d-good-initial -points gives us our 4 sets of d initial 
values on the circle of radius 3/2. The sets have the property that elements of the 
same set are mapped to points with approximately the same argument, and elements of 
distinct sets are mapped to points in different quadrants of the target space. There are 
several different ways to accomplish this, and other methods can slightly decrease the 



number of operations required in iterate-PLM. Please refer to the remarks in Section 
for more details. 

function choose-4d-good-initial -point s(/) 

Let N = 676 deg/. 
For A; = 1 to iV 

Let cvk = fe^-'^/^. 
End For. 
For j = 1 to 4 
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Let Zj be the union of the cja: for which wigfiuo]^ < 27r/j and arg f{ujk+i) > 

2vr/i. 
End For. 
Return(2). 



From our sets of initial points Zj, we obtain our approximate zeros via the routine 
iterate-PLM. This applies the path lifting method to / for an appropriate number 
of steps. For simplicity, we present a scalar version here. However, there are well- 
known methods for evaluating the same polynomial / at m different points (m > d) in 
O (rn{log d)^^ arithmetic operations; refer to section 4.5 of |[BM75|| . When computing 
the complexity in Section P?7| , we assume such methods are used. Note that this algo- 
rithm can be easily implemented on either a vector or parallel computer, so that one 
can iterate the d elements of Zj simultaneously. Also note that if Smale's mean value 
conjecture holds, the extra iteration to obtain z is not needed; zn is good enough. See 
the remarks following Lemma O. 



function scalar-iterate-PLM(/, ZQ,j, r) 



Let wo = \fizo)\e^^'^^. 
Let h = 1/27. 

LetiV= ^§4^ 
Llog2 (26/27) 
For n = 1 to iV 

Let Wn = — h)Wn-l- 

f {Zn-l) - Wi 



Let Zj 
End for. 
Let z = zj\f 
Return(i). 



Zn-l 



f'{Zn-l) 

f{zN) - rei^'l^ 

f'izN) ■ 



The next routine takes the output of iterate-PLM and uses the a function (defined 
in Lemma p..3|) to remove those elements which are not approximate zeros. Although 
the test a < 1/8 is sufficient, it is not a necessary condition. However, if the image of 
the initial points Zj lie in a "good quadrant" , we know by Lemma p.5| below that we 
will have a < 1/8 for at least half of them, and so we are not in danger of discarding 
too many points. 



function select-approx-zeros(^/', 3^) 

Let ;f = 0. 

For j = 1 to #{y) 

If a^ivj) < 1/8, then let X = X U {yj}. 
End For. 
Return(;f). 
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Once we have found approximate zeros for at least half of the roots of fi (by applying 
iterate-PLM to at most 4 sets Zj), we further refine them by "polishing" with regular 
Newton's method. As in iterate-PLM, fast polynomial evaluation techniques can be 
used, but we present a scalar version here for simplicity. 



function scalar-polish-roots('0, Xo , r) 

log2log2 (64^(7/4)'=') /i 



Let M = 1 + 
Forn = 1 to M 

Let Xji = Xji—i I N • 
End For. 



Return (xm)- 

Finally, we remove from the set of approximations that we have found any points 
which approximate the same root of ■0. We do this by taking the approximation which 
minimizes Iz/^l, and then making a pass through the rest of them and accepting only 
those which approximate different roots from the ones previously accepted. 



function weed-out-duplicates('0, W) 

Sort W so that \tp{wi)\ < \tp{w2)\ <■■ ■ < \tp{wn)\- 
Let V = {wi}. 
For J = 2 to #(>V) 

If \wj -v\> 3\'tp{wj)\/\'tp'{wj)\ for all veV, then let V = V U {wj}. 
End For. 
Return(V). 



At this point we have found approximations to at least half of the roots of fi. We 
divide them out to obtain a new polynomial /i+i of smaller degree, using a standard 
technique involving the finite Fourier matrix. 



function def lateC^^, V) 

Let n = deg — #V. 

Let Lo = e2'^V('^+i). 
#v 

Let p{z) = Yl{z-Vk). 

k=l 

Let M"^ be the inverse of the n x n Fourier matrix. That is, mj ^ = ui^^^ /n. 
Let V be the column vector with entries p-j = — ; — r^, 7 = 0, . . . , n. 

Let Q = V. 

Letqiz)=j:]=oQjZ^- 
Return(g). 
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We now state our main theorem, which essentially says that the algorithm just 
presented works: 

Theorem 2.1. Given a monic polynomial (j){z) of degree d andt > 0, the algorithm 
presented in this section will always terminate with d points Ai, . . . , which satisfy 

Uiz) - {[{z - XM < e. 

For (f) G Vd{^), the worst case arithmetic complexity of this algorithm is 

0{d{\ogdf\ logel + d\logdf). 



Corollary 2.2. Let (f) and e be as in the theorem, and denote the roots of (j) 
by C,j, represented with multiplicity. Then the algorithm can be used to produce the 
approximations Xj satisfying 

with a worst- case arithmetic complexity of 

0[d\logdY\loge\+d\logdf). 



Proof. This is an immediate consequence of the fact that if / and g are in 'Pd(l) 
with 11/ — g\\ < {e/8dY, then the roots of / and g are at most e apart (See ||Kim89b|| ). 
□ 

We prove the theorem as a series of lemmas in the following sections. There is a 
rough correspondence between the sections and the routines in the algorithm. Finally, 



we summarize all of these lemmas and give the proof in Section 2.8 



2.2. Selection of initial points. For each intermediate polynomial fi, we need 
to select four sets of (deg fi) points at which to begin our iteration. These must be 
chosen so that the elements of each set map very near the same point in the target 
space, and that the images of elements in successive groups are approximately —turn 
apart. This can be accomplished either by evaluating fi at a large number of points 
spaced around a circle in the source space, and then selecting from those, or by taking 
a much smaller number (perhaps as few as Ad) of points and adjusting them with either 
a standard or modified Newton's method. The arithmetic complexity of either comes 
out much the same; we opt for the former method because of its conceptual simplicity. 

The following lemma gives us bounds on how much the argument in the target 
space can vary between points around a circle containing all the roots in the source 
space. Using this, we see how many points are required to obtain our "good" points. 

Lemma 2.3. Suppose all of the roots C of a polynomial f of degree d lie in Dji{0), 
and let 
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Then 



arg 



m+lj 



< 



n 



Proof. This argument appears in | Ren87|| (Lemma 7.1), although in a somewhat 
different form. We present an adapted version here. The idea is quite simple: since 



arg- 



m+l) 



arg ■ 



n- 



E^m+l ~ Ci 
arg — , 

.1—1 Si 



we merely bound the angles in the source space and add them up. 




Fig. 2.1. Calculation of the upper bound on the variation of argument in the target space. This 
picture is in the source space. 



Note that each root \Ci\ < R and \uj. 



ml ^ 2R, so 



arg ■ 



^m+l — Ci 



Ci 



< 2 arctan 



2/?cos^ - R 

na 



< 2 arctan ■ 



2i?sin^ 

na 



2i?(cos 



< 



An 
nd 



Refer to Figure 2.1. Summing the d terms gives the result 

for h 



From Theorem 1.5 



1/27 we require that our initial points be within 
sin~^(/i/2) < 7r/169 of the central ray, so we start with 676c/ points equally spaced 
around the circle of radius 3/2. We then evaluate the polynomial at each of them, and 
select four sets of d points whose arguments are closest to 0, |, tt, and respectively. 
For each of these, we take the initial target point wq to be the projection of its image 
onto the real or imaginary axis. 

Remark. Some amount of computation can be saved if we use the same wq for 
all d elements of a given set of initial points, rather than just points with the same 
argument. This would make the computation of the tf„ a scalar rather than a vector 
operation, that is, Wn would only need to be computed once for each group of d points. 
However, in order to do this, we must ensure that the images of each Zj in the same set 
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have approximately the same modulus as well as argument. This is best accomplished 
using some sort of Newton's method. 

If an initial Newton's method is used, one can also choose a much smaller number 
of trial points cuj. Such a method should converge well, since all the critical points of 
f{z) — wq are inside the Br, while the roots of f{z) — wq and the Ui are well outside. 



2.3. Iteration of the path lifting method. In this section, we analyze the 
behavior of applying the path lifting method to a single well-chosen initial value zq. 



We first prove the theorem as promised in Section |1.3| , and then we show that after 
a specified number of iterations, the result will be an approximate zero. Before doing 
this, we shall state the relevant special case of Theorem 3.2 of ||Kim88|| , which gives a 
lower bound on how far a point moves under Newton's method. 



Lemma 2.4. Let z 



g{z)/g'{z), with r = \g{z)\/Rg{z) < 0.148, where Rg 



is the radius of convergence of g^ as a power series based at z. Then 



< B(r) 



where 



B(r) 



2r- 



Given this lemma, the proof of Theorem |1.5| is quite simple. Note that each step of 
the iteration in Theorem corresponds to Newton's method applied to f{z) — Wn+i 

at Zn—l- 

Theorem |1.5| . Suppose that the branch of the inverse f^J is analytic on a wedge 
y^A,woj ''^'ith < A < Tx 12, and let h < Suppose also that \ f{zo) — Wq] < h\wo\/2, 

and define 



Wr. 



h^wo, 



f{Zn) - Wn+1 
f'i^n) 



Then |/(2„) - w„| < h\wn\/2 and z^+i E /^,/ (^a,wJ ■ 

Proof. We shall prove this by induction. All that is required is to establish the 
conclusion, given that — Wn-i\ < h\wn-i\/2. 

Note that 



- Wn\ < \ f{Zn-l) - Wn-l\ + \Wn-l - Wn\ < 3h\Wn-l\/2, 

and that 



Rf-Wn{Zn-l) = Rf{Zn-l) > \Wn\smA- h\Wn\/2. 



Since h < ^^fa-, we can apply Lemma |2]J with g{z) = f{z) — Wn and r = 3/37 to obtain 



f{Zn) - Wn 




f{Zn) - Wn 




f{Zn-l) 


- Wn 


Wn 




f{Zn-l) - Wn 




Wn-l(l 


-h) 



< B(r] 



3h 
2(1 - h)' 



Because 5(3/37) < (1 — h)/3, we have our conclusion. □ 
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Remark. The denominator of 19 in the upper bound on h can be relaxed only 
slightly in this proof. If we take h < f^, then k(A) is a monotonically increasing 
function, with 18.3096 < k{A) < 18.895. Also, if A > 7r/2, we can take h = 1/19. 

We apply the path lifting method until we have an approximate zero, at which time 
we can use the standard Newton iteration. Since we are interested in finding r-roots of 
/, we stop the iteration once we have an approximate zero for / — Tj^- Note that, as 



discussed in Section 2A, this translation is necessary, since there are no approximate 
zeros in a neighborhood of a multiple root. Having approximate zeros for the perturbed 
polynomial gives us the following bound on the number of iterations required, and 
enables the weeding out of duplicates. 

Lemma 2.5. Let h < with Wi and Zi as in Theorem \1.^. If 
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N 



log2(^/ko|) 
l0g2(l-/l) 



77 \ tfl/GTh 



then Zn is an approximate zero for f — Tj^. Furthermore, if z = zn 
we have a{z) < 1/8. 

Proof. For notational convenience, set r = Tj^- 

By assumption, f~^ is analytic on the wedge yVA,wo, and the initial point zq is close 



enough to wq that we can apply Theorem \L.5[ Note that if is as specified, we have 

r 



since |wAr| = (1 — h) \wo\ by definition. Also, Ifiz^) — wn\ < h\wN\/'2. Refer to 



Figure 2.2 



f(z^) in here 




Fig. 2.2. The locations of wn and f{z]sf) in the target space when zq is in a "good quadrant". 



In order to conclude that zn is an approximate zero, we need to show that we 
have \f{zN) — t\ /Rf{ziy) < 1/10, where Rf^z^) is the radius of convergence of f~\ 
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Consider /{zn) in some circle of radius h\w\/2, where argw = argwo, and \w\ > r. 
Then we have 

\f{zN)-r\^ \w\ + h\w\/2-T _ 2 + h 



\w 


\ + h\w\/2 


- T 


\w\ 


sin A — h 


w\/2 



Rf{zN) ~ \w\smA-h\w\/2 2smA-h \w\{smA- h/2)' 



The maximum of this quantity occurs when \ w\ is as large as possible, which in our case 
is \w\ = r/(l — h). Hence 

I fizn) — t\ 3h 

'^„7 , ' < 7-—- r < 3/37 < 1/10. 



By Lemma |1.4] , zj^j is an approximate zero of f{z) — r. 

In order to ensure that a < 1/8, we need to apply Newton's method once. Set 

f{zN) - T 

z — Zn — ^~• 

Then this point satisfies 



\f{zN) 

and so 



< fi(3/37) < 0.31271, 



\f{z)-T\ < 0.31271 {\f{zN)-WN\ + \wN\-r) < 0.31271-^^ < ^^'""^ 



2{l-h) - 39 



This gives us 



\f{z) -t\ ^ T sin A/39 ^ ^^gg^ 



Rf{z) rsin A - rsiny4/39 
As a consequence of Lemma |L^, we have a{z) < 1/8. □ 

2.4. Refinement of the root approximations. The routine iterate-PLM out- 
puts a set of points y = {yi,y2, ■ ■ ■ ,yd} which may be approximate zeros for ip{z) = 
f{z)—f. We can use (see Lemma PT^) to choose those yi that are indeed approximate 
zeros; we discard those yi for which a^{z) > 1/8. As a consequence of Lemma |1]^ and 
Lemma those yt which started in the "good sector" will not be discarded. 

In addition, we also want to ensure that we approximate each root only once 
(counted with multiplicity). Lemma gives us conditions which allow us to weed 
out any duplicates; the proof relies on a variant of the Koebe Distortion Theorem 
which we quote here from [|Kim88|] , Lemma 3.3. 

Lemma 2.6. Let g be univalent on Dji{z). Then for any s < 1, we have 

BtR {g{z)) C g (BsRiz)) C BuR (giz)) , 
where t = s\g'{z)\/{l + s)^ and u = s\g'{z)\/{l - s)^. 
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Lemma 2.7. Suppose thatyi andy2 are approximate zeros ofip{z), with \%p{yi)\ > 
\4'{y2)\, and R^iijj) > Vlipiiij), where R^iyj) is the radius of convergence of ip'K Then 
yi and y2 approximate the same simple root C, of ip if and only if 



Proof If yi and 1/2 approximate the same simple root, then we have ipy_^ = ip^^. 
Since - ip{y2)\ < 2|V'(|/i)|, we have ?/'(|/2) e ID)2|V'foi)|(V'(l/i))- Thus we can apply 

Lemma |2.6| with g = which is univalent on the disk of radius R = 12|?/'(?/i)|. 

Taking s = 1/6, we have ip"'^ (D2jV'(?yi)|('^(yi))) contained in the disk of radius 
about yi, so \yi - 7/2! < i^^^- 



25|V'fei)| 

For the other direction, we apply the Koebe ^-Theorem (or Lemma p.6| with s = 



1; it is the same thing) to see that ipy^ (®i2|V'(?/i)|(V^(l/i))) contains the disk of radius 
3\ip{yi)\/\4''{yi)\ about yi. Thus, if the distance between yi and y2 is less than this 
amount, we must have ip~_^{ip{y2)) = y2, that is, yi and y2 approximate the same root 
of tp. □ 

Remark. Note that < 1/8 is not sufficient to imply that the hypotheses of 

Lemma are satisfied, since this only gives R^{z) > 4:\ip{z)\/3. However, if Newton's 
method is applied to such a point at least 3 times, the value of\ip\ will decrease by at least 
1/128, and so the lemma can be applied. Since we need to "polish" the approximations 
with Newton's method in order to control the error in the deflation, we do that before 
weeding out the duplicates. The total number of iterations of Newton's method required 
is calculated in Section ^]6|, but it is greater than 3 in all cases. In practice, one should 
probably perform the weeding in the routine select-approx-zeros, in order to avoid 
polishing points which will be discarded later. 

2.5. Deflation of intermediate polynomials. Here we compute an explicit 
bound on the error introduced by the deflation step. We start with a polynomial ip 
which has roots {C,j}j=i,...,d, and a set of approximations to these roots which we de- 
note by {vk}k=i,...,n with n < d. We then use polynomial interpolation via the discrete 
Fourier matrix to obtain a polynomial q of degree d — n so that 

n 

ip{z) ^ p{z)q{z), where p{z) = Y\_{z — Vk)- 

k=l 

In this section, we estimate Wip^z) — p{z)q{z)\\, as well as the accumulated error in 
repeating this process until q{z) is a constant. Recall that the finite Fourier matrix Ai 

has as its j, k^^ entry the jk^^ power of a primitive d + l^*-root of unity. That is, 

m,- fc = u^'' = e^^'^^/^^+^) j,k = 0,...,d 

Then one can readily see that if f{z) = Y.'j=o djZ^ and A is the column vector of the 
coefficients of /, then the product M.A will be the vector B whose j^^ entry is the 
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value of f{oJ^). Also, given the values of / evaluated at the powers of a;, we can easily 
compute the coefficients of / as the product M.^^B. 

Our first lemma gives us estimates on the size of the error caused by a single 
defiation. As is common with this sort of thing, the proof is neither terribly entertaining 
or enlightening. 

Lemma 2.8. Suppose iIj{z) = 11^=1 (-2 - ^j), with < 3/4, factors as ■iIj{z) = 
P{z)Q{z), where deg(P) — n < d and deg(Q) — m — d — n. Let p{z) — n"=i(-2 — vj), 
with l^j — Vj\ < S < and define q{z) and r{z) by iplz) — p{z)q{z) + r{z), where q is 
found by polynomial interpolation as described above. Then 

\\P{z) - p{z)\\ < 8n(5(7/4)" 
\\Q{z) - q{z)\\ < 8m5(7/4r 
||r(^)||< 8d5(7/4)^ 



Proof. Let cu = 6^^^'/*^"+^^ and let B be the vector with entries bj 
for < j < n. Note that 



Also note that 



k=l k=l 



n(^^-c.) 1-n 



yk=i 

/ n 



= n^-c^) 1-n 1+ 



k=l 



- Vk 

- ik 



4<k-e.|<4 



and 



ik - Vk 



< 45. 



Since ^5 < < we have (1 + 45)" < 8nS + 1, and so 

\bj\ < (7/4)"((l + 45)" - l)< 8n5 (7/4)" . 

Thus 

\\P{z) -p{z)\\ = \\M-^B\\ < \\B\\ < 8nS{7/4y 



We have a similar computation for the bound on HQ — Let rj — e^'^'/("*+^) and 
let C be the vector with entries Cj — Q{r]^) — q{r]^), j — 0, . . . ,m. Then 



p(r^j) p(r^i) p(r^i) \^ p(,^i) J y ^J-^ rj^ - 
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As in the case of P{z) —p{z), we obtain 

\\Q{z) - qiz)\\ = \\M-'C\\ < \\C\\ < 8m6 (7/4)" . 

Finally, for the bound on ||r(2;)||, note that r{z) = P{z)Q{z) —p{z)q{z). If we write 
= Yfj=Q^jZ^-, P{z) = Yfj=oPjZ^ I and so on, then we have 



Pj-kQk — '^Pj-klk = '^ylk{Pj-k — Pj-k) — Pj-kiQk — Qk))- 



k=0 k=0 k=0 

Since < 3/4, we have the following crude bounds on the coefficients of P and q: 

\Pj-k\<i7/Ar and |g,|<(7/4r. 
Combining this with the bounds on ||P — p|| and HQ — gH, we get 

\rj\ < ^ (8m5 (7/4)" (7/4)" + 8n5 (7/4)" (7/4r) 

k=0 

< 8d^6 (7/4)*" . 

□ 

Now that we have bounds on the error in one step of deflation, we can bound the 
error introduced by repeated deflation. We assume that our initial polynomial / has 
roots in so that we can ensure that the roots of the subsequent polynomials fk 
remain in D3/4 as required by Lemma |2.8| . 

Lemma 2.9. Suppose f{z) = IYj=i{z - Q), with < 1/2. Let /o = /, and define 

fkiz) = Pk+iiz)fk+iiz) + rk+iiz) 

for < k < m — 1, where pj, fj, and rj are determined by polynomial interpolation 
as in Lemma |^7^ , with fm,{z) = 1. Suppose also that degr^+i < deg = dk, and that 
\\rk\\ < min (^^^ ^{4/7^+^) . Then 

\\f -P1P2 ■ ■■Pm\\ < /i- 

Furthermore, if we have dk < dk-i/2 and fj, < (7/16)'^, we need only require that 
\\rk\\ < I2{A/7Y+^ 

Proof Note that / - piP2 ■■■Pm=Pi- ■ ■Pm'^m-i + • • • + Pir2 + n, and so 

\\f -P1P2 ■ ■■Pm\\ < \\Pl ■ ■ ■ Pm-l\\\\rm-l\\ + • • • + |bl|||k2|| + 

First, we show that ■ • - pfeH < (7/4)<iegpi---Pfc^ shall use induction to show that the 
roots pi - ■ -pk are always in the circle of radius \ + < 3/4. On the circle \z\ = | + ^, 



we have 



deg rfe /7 1 

\fk^^{z)-pk{z)fk{z)\ < \\rk\\ E < ,,.,,0^ M ^ ^ 1/^-1(^)1, 



i=o 



{AdY{2d -k) - {Ad) 
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and so by Rouche's Theorem the roots of Pkfk he inside Thus the roots of 

P1P2 ■ ■ -Pk also he in this disk, and so the coefficients are less than (7/4)'^°spi '^'=. 
Since degpiP2 ■ ■ ■ Pk < d — {m — k) , we can conclude that 



m— 1 
k=l 



■■Pk\\\\rk+i\ 



< 

< 
< 



m— 1 
k=l 

i(T/4) 

/i. 



d+2 



(4/7)''+'// 



If we require that dk < d/2^, that is, we find at least half of the roots at each step, 
then we can relax the restriction on ||rfc|| somewhat. If we have ||rfc|| < then we 
can use Rouche's theorem as before to show that the roots of Pkfk are in the disk of 
radius Ck 



1,1 v^fc 
2 "T 2 ^n=l 



2-2' < 3/4. Note that for |z| 

1 



\fk-iiz)\ > (Ck-Ck-iY 



> 



.2-22'= 

Applying almost the same argument as before, we have 



Ck, we have 
> 1/4° 



\fk~i{z) -Pk{z)fk{z)\ < \\rk\\ J2 



j=0 



\z\' <^,<\fk-i{z)\. 



Since ||r,t|| < /i(4/7)^ < j^, we are done 



□ 



2.6. Controlling the error. In this section, we compute the size of r that we 
can use to ensure we have an e-factorization of tp. The following very simple proposition 
shows that if the norms of two polynomials are close, so are the norms of the rescaled 



versions. This gives us the relationship between e and the number fi used in Lemma 2.9 



Proposition 2.10. Suppose that cj) is a monic polynomial with all its roots in Hr, 
so that the roots of f{z) = {2R)~'^(j){2Rz) are in D1/2. Then if\\f — p|| < e, we have 



{2RYp 



\2R 



< {2RYe. 



Proof. Let f{z) = J2(^jZ^ and p{z) = J^bjZ^. Then (f){z) = J^i'^RY ^Ojz^ and 
{2RYp (2^) = j:{2RY-^bjzt Since max | E(2i?)''~^(aj - bj)\ < {2RY\aj - bj\, we have 
our claim. □ 

Our input polynomial (p is in 'Pcz(l) and / is be the rescaled polynomial as in the 
previous proposition, so i? = 2. Then an e-factorization of corresponds to an e/4'^- 
factorization of /, so we take /i = e/4'^. 

In order to properly approximate the roots of /, we need to ensure that the remain- 
der (as in Lemma |2]9| ) at the k^^ step satisfies ||rfc|| < 2r, where r = {4/7Y~^^^/2 = 
32e/7'^+3. 
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At each stage, we translate fk by r, and ensure that the error introduced by the 
deflation of the translated polynomial is no more than r. By Lemma p.8| , we need the 
root distance between translated polynomial and the deflated polynomial to satisfy 

- 8d \7 
Then we will have 

Ikfcll = ll/fe-i -PkfkW < Wifk-i - r) -pkfkW + Wifk-i -r)- fk-i\\ <T + T. 

In order to achieve the root distance less than 6, we apply Newton's method to 
the approximate zeros found by the routine iterate-PLM (see Lemma |2.5| ). Since each 
point z is an approximate zero to the root ^ of the translated polynomial, we have by 



Deflnition L2 



5<8(- k-el- 



.2. 

Thus, iterating Newton's method log2 log2(8/(5) times, as is done in polish-roots, will 
give the desired result. 

2.7. Arithmetic complexity. In this section we count the number of arithmetic 
operations involved in using the algorithm to obtain an e-factorization of a polynomial 
in Vd{l). 

In the main routine, we rescale the polynomial and then invoke get-half -roots- 
and-def late at most log2 ci times, since at least half of the roots are found in each 
call. 

• The cost of rescaling is 2d multiplications. 

• Each call to get-half-roots-and-def late calls choose-4d-good-initial-points 
once, and makes at most 4 calls to each of iterate-PLM, select-approx-zeros, 
polish-roots, and weed-out-duplicates, and one call to deflate. As before, we 
denote the degree of the input polynomial by d and the degree of the A;*^ intermediate 
polynomial by dk. 

o The routine choose-4d-good-initial -points involves evaluation of fk at 676dk 
points, and 676dk comparisons. Since we can evaluate fk at m points with a 
cost of O (m{\og dk)"^^ operations (see ||BM75|| ), this gives a total of 0(^dk(i.ogd 



i-kf 



operations. 



o For iterate-PLM, each iteration evaluates Zn — '^''^//(z"')"^^ dk points, which costs 
0(^dk{logdky^- This is done times, where A^ < 271og(^). Note that since 
the roots of / are in D3/4 and \zq\ = 3/2, we have \wo\ = \ fk{zo)\ < (9/4)'^*. Since 
r = 32e/7'^"'''^, we have A^ = Ci{d+\ loge|) for some constant Ci. This gives a total 
of 

C(((i+|loge|)(4(log4)^ 

operations. 
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o For select-approx-zeros, we need to evaluate a at points, which requires 
evaluation of all the derivatives of /. This requires 0(^dl{\ogdk)'^^ operations. 

o The routine polish-roots performs M iterations of Newton's method at dk points. 
Since M = C2(}ogd + log | loge|) for some constant C2, we have a total operation 
count of 

C»((4(log4)')(logrf + log|loge|)). 

o Weed-out-duplicates requires a sort of at most 4 points, which costs 0{dk log d^), 
and evaluation of / and /' at dk points. Thus the total cost of this routine is 

o(4(log4)'). 

o Lastly, deflate costs (9 (^4 (log 4)^) operations. 

Thus, the overall cost of each call to get-half-roots-and-def late is dominated by 
that of iterate-PLM, and is at most 0(^4(log4)^('^ + I loge|)j operations. Since 
4+1 < 4/2, we have a total cost of at most 

0[d\logdy + d{\ogd)^\loge\) 

operations to obtain an e-factorization of the polynomial. 

2.8. Summary and proof of main theorem. At this point, we have actually 
already proven Theorem but we would like to tie together the various steps involved. 
Just to refresh your memory, this theorem says, in essence, that our algorithm always 
produces an e-factorization with the stated complexity. 

Recall that the algorithm performs the approximate factorization in stages using 
the routine get-half-roots-and-def late; at the fc*^ step, we produce a function 
and sets of approximations so that 

AiSA 

where A = Uj<fc ^j- I'^ order to prove Theorem ^.1| , we need to show that the number of 
approximations found at the k^^ stage (#Afc) is at least (deg/j!c_i)/2, and that \\f{z) — 
fk{z) Y[xi£A{z ~ < 6/4*^. Given this, the complexity calculation in the previous 
section will apply, and we shall have the result. 

First, note that as a consequence of Lemma |1.6| , there will always be a quarter- 
plane in the target space on which there at least dk/2 branches of f^^ are defined. 
This means that if we start with 4 points zj which are well-spaced in the source 
space (so that each sheet in the target space is represented), then for at least half 
of them there will be a branch of the the inverse which is defined on the entire 

j 

quadrant, and 7^ f~\ For these points zj, if we ensure that f{zj) is close enough 
to the center line of the quadrant, at least half of them will satisfy the hypothesis of 



Theorem |1.5| and so the routine iterate-PLM will produce approximations to each of 
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the corresponding dk/2 roots of fk (with multiphcity). Such initial points Zj will be 
produced by choose-4d-good-initial -points, as was shown in Section |2.2| . As a 



consequence of Lemma p.5| , the good approximations iterate-PLM are approximate 
zeros of ip = fk — T, with < 1/8. As was discussed in Section application of 
the routines select-approx-zeros and weed-out-duplicates will select exactly one 
representative for each approximated root of tp, giving at least dk/2 such approximations 
Xj. This selection is necessary since some of the initial Zj which did not satisfy the 



hypothesis of Theorem |1.5| may still have converged. 

As was shown in Section |2.6| , the approximations produced yield an e/4'^-factor- 
ization, since the each the Aj are made sufficiently close to the roots of ip by the routine 
polish-roots, and — fk\\ is sufficiently small. This completes the proof of Theo- 
rem ^.1| . 
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